A numerical technique for studying topological effects on the thermal properties of 

knotted polymer rings 
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The topological effects on the thermal properties of several knot configurations are investigated 
using Monte Carlo simulations. In order to check if the topology of the knots is preserved during the 
thermal fluctuations we propose a method that allows very fast calculations and can be easily applied 
to arbitrarily complex knots. As an application, the specific energy and heat capacity of the trefoil, 
the figure-eight and the 8i knots are calculated at different temperatures and for different lengths. 
Short-range repulsive interactions between the monomers are assumed. The knots configurations are 
generated on a three-dimensional cubic lattice and sampled by means of the Wang-Landau algorithm 
and of the pivot method. The obtained results show that the topological effects play a key role for 
short-length polymers. Three temperature regimes of the growth rate of the internal energy of the 
system are distinguished. 

PACS numbers: 



I. INTRODUCTION 

Topological polymers are a fascinating subject by 
themselves. Moreover, the interest in this topic is mo- 
tivated by possible industrial applications, for example 
because knots in polymer materials are able to influence 
the viscoelastic properties [1], and due to their relevance 
in biology and biochemistry, see for instance [2l-llQ| . Poly- 
mer rings may be found in nature or in artificial materials 
in the form of knots that can be additionally linked to- 
gether. In this work we explore the case of the statistical 
mechanics of a single polymer knot. Various topological 
configurations are taken into account: The unknot, de- 
noted according to the Alexander-Briggs notation Oi, the 
trefoil 3i, the figure-eight 4i and the knot 8i. In partic- 
ular, we study the thermal properties of the above knots 
by computing numerically their internal energy and heat 
capacities. The most important problem of numerical 
simulations of polymer knots is to preserve their topol- 
ogy while their configurations are fluctuating randomly. 
To this purpose, several methods have been developed 

una. 

Most of them, are based on the Alexander or 
Jones polynomials, which are rather powerful topologi- 
cal invariants able to distinguish with a very high degree 
of accuracy the different topological configurations. The 
main drawback of these polynomials is that their calcu- 
lation is very expensive and time consuming from the 
computational point of view. 

Here we propose a new strategy that consists in the 
following. First, we apply the pivot algorithm to change 
a given knot configuration. Starting from a seed knot, for 
instance those given in ref. [18| , we change it at randomly 
chosen points by applying the pivot algorithm [l9| . After 
each pivot move, the difference between the old and new 



configurations forms a closed loop. Around it, an arbi- 
trary surface is stretched, whose boundary is the closed 
loop itself. The criterion to reject changes that destroy 
the topology of the knot is the presence or not of lines 
of the old knot that cross such surface. If these lines 
are present, the trial pivot move is rejected. This com- 
bination of pivot algorithm and excluded area (PAEA) 
provides an efficient and fast way to preserve the topol- 
ogy that can be applied to any knot configuration, inde- 
pendently of its complexity. For pivot moves involving 
a small number of segments the method becomes exact. 
This technique may be employed in the study of the ther- 
mal and mechanical properties of polymer knots. Here 
we have limited ourselves to the computation of the in- 
ternal energy and heat capacity in the case in which the 
monomer interactions are of very short range, but the 
extension to any other interactions or to the calculations 
of more complicated physical quantities is possible. The 
configurations of the unknot, trefoil, figure-eight and 8i 
knot are modeled on a cubic lattice. The sampling of 
the canonical ensemble is achieved by using the Wang- 
Landau algorithm [20] at different temperatures. 

The rest of this paper is organized as follows. In 
Sec. II, we describe the sampling and simulation meth- 
ods of the polymer knots. Our strategy to distinguish 
possible violations of the topology of a given knot after 
a random change of its trajectory is explained. The 
calculation of the specific energy and the specific heat 
capacity for various polymer topologies and for different 
polymer lengths is presented in section III. In section IV 
we draw our conclusions and discuss the open problems 
together with new possible directions of our research. 



II. MODEL AND COMPUTATIONAL METHOD 
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In the current work, we will study the topological ef- 
fects on the thermal properties of polymer rings. Four 
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FIG. 1: In this work different knot types will be studied, 
namely the unknot, the trefoil, the figure-eight and the 81 
knot. 



types of polymers will be considered on a 3-dimensional 
cubic lattice, namely the unknot U, the trefoil T, the 
figure-eight T and the knot 81, see Fig. [TJ A general 
knot will be denoted with the symbol JC. To simulate the 
trivial polymer ring we generate a loop starting from 
a random point and performing a self-avoiding walk on 
the cubic lattice. The walk is stopped when it crosses 
itself, forming a loop. We generate many loops of this 
kind and we retain only those with a given length L. 
For nontrivial topologies, we choose a different strategy. 
We begin from a seed conformation and then we act re- 
peatedly on it with the pivot method [l9| to get random 
conformations. The elements of the loop that will be 
changed by the pivot move have a fixed number of seg- 
ments TV = 4 or TV = 5. The elements begin at a random 
point in the knot, let's say Nq and end up at the point 
number TVo + 4, or TVo + 5. TVo can range from to L — TV 
(L is the length of the polymer). We check the topol- 
ogy of the configuration after every pivot move, to be 
sure that the move will not change the topology of the 
given knot. Usually, this check is achieved by projecting 
the knot along some direction on a two dimensional plane 
and the relative Alexander polynomial is computed. This 
method has been pioneered in refs. Here we choose 
another strategy, namely on a given knot JC we perform 
a pivot move. In that way a new knot JC' is obtained. A 
schematic picture of the situation can be found in fig. [2l 
In order to verify that the topology of the knot has not 
been changed by the pivot move, we check that no lines of 
JC cross the area spanned by the small loop AJC = JC' — JC 
including its border. Otherwise, this could imply that a 
crossing of the lines of JC has happened. The definition 
of the area of AJC in a numerical simulation is a difficult 
problem of computer graphics. However, if the number 
of segments involved in the pivot moves is small, it is pos- 
sible to classify all the geometries of AJC and to stretch 
a suitable surface Xak: around AJC. Next, these data, 
namely all geometries of AJC and the coordinates of the 
relative surfaces Eak;, are supplied to the code, which will 
then be able to detect the occurrence of possible lines of 
JC crossing Eak: or the part of its border AJC — JC f] A/C, 
see fig. [2] for a schematic representation of what these 
symbols mean. If these lines are present, the trial pivot 
move is rejected. Let us note that the strategy described 
above could be considered as an implementation of the 
topological constraints based on the Gauss linking num- 
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FIG. 2: This figure represents in a schematic way the meaning 
of the symbols used to denote the element JC f] AJC of AJC 
shared with the knot JC and the element AJC — JC f] AJC of 
AJC which is not shared with JC. 

ber modified to cope with the present situation, in which 
the Gauss linking number becomes not suitable because 
JC and AJC share part of their trajectories. While the 
Gauss linking number is a weak topological number, our 
method turns out to be very efficient when applied to 
pivot modes with segment lengths TV = 4 and TV = 5. In 
both cases the geometry of the loop AJC is simple enough, 
that it is possible to analyze all situations in which the 
crossings of lines could violate the topology of JC and 
to provide an exact criterion for preserving this topology 
without increasing exceedingly the simulation time. Only 
when TV > 5 it may occur that the trajectory of the old 
knot crosses the area Xak: an even number of times. In 
this case, due to limitations similar to those of the Gauss 
linking number, allowed moves can be rejected. If this 
happens very frequently, the computational time may in- 
crease. However, this is not a serious drawback because 
the above procedure of detecting topology changes is ex- 
tremely time efficient. The necessary cpu time depends 
essentially only on the number of segments L, indepen- 
dently on the complexity of the knot. In the following, 
we will restrict ourselves to the cases TV = 4 and TV = 5. 
In order to calculate the internal energy of the polymer 
per unit of length (E(f3))/L, where f3 = T _1 denotes the 
Boltzmann factor with ks = 1 (natural thermodynamic 
units) and the heat capacity 

C(T) = ^{{E{P?) ~ (W)» 2 ) (1) 

the density of states of polymers is needed. To that pur- 
pose, the Wang-Landau (WL) algorithm is used for sam- 
pling the polymer conformations. In the case of the tre- 
foil, figure-eight and 81 knot, before applying the WL 
algorithm, the system has been equilibrated, since we 
start from a given seed which is not in a random confor- 
mation. The fluctuations of the radius of gyration Rq 
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FIG. 3: M is the number of sites of the trefoil in which the 
direction of the chain is not changed. 



are employed to judge the degree of equilibration. These 
fluctuations are computed according to the formula (2l| : 



AR%(t) 



\ 



(2) 



Here t is the sweep time, {R%)t is the average radius 
of gyration at time t. The equilibrating procedure is 
stopped after the time dependent fluctuation of the ra- 
dius of gyration (ARq) oscillates around a constant. To 
crosscheck the equilibration of the nontrivial knot con- 
figurations, we use also another method, that consists in 
counting the number of sites in which the direction of 
the chain is not changed [22j. This number is a random 
variable, but at equilibrium its variance around the aver- 
age is small. The results obtained with both techniques, 
that based on the fluctuation of the gyration radius and 
that based on the changes of direction of the segments, 
are in agreement. Fig. [3j which is related to the latter 
method, shows that the system equilibrates after a few 
tens of thousands of pivot moves in the case of the trefoil 
T. Similar considerations are valid for the other knots 
considered in this work. In total the system is subjected 
to one million of pivot moves before treating it with the 
WL algorithm in order to compute the energy (E(/3)) 
and the heat capacity C(T). 

Now we briefly introduce the WL algorithm. This 
method is considered as a self-adjusting procedure for 
obtaining the density of states. In the present context 
the states are distinguished by the number i of closest 
contacts between the polymer monomers, where i takes 
positive integer values. Following a previous work [23| , 
all initial values of the density of states are firstly set to 
be equal to one, i. e. = 1 for all the states. Then the 
pivot method is used to visit the next trial state i' . The 
energy of the trial state will be denoted by era, where m 
is the number of the closest contacts of the monomers, 
while e is the interaction energy of each monomer pair. 
e > and £ < correspond to repulsive and attractive 
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FIG. 4: The density of states is obtained by WL algorithm. 
The results for the unknot (triangles) , the trefoil (circles) and 
the figure-eight (rectangles) are represented. 



case respectively. We assume in the following that e > 0. 
The probability of translating from the state i to %' is 



p(i — i') = mm[l, 



(3) 



Once an energy state is visited, the corresponding den- 
sity of states of this state will be updated by multiply- 
ing a modification factor [20]. After the n— th sweep, 
the modification factor fa is reduced according to the 
equation fa+i = \ffa- The above process is repeated un- 
til the modification factor reaches the predefined limit 
f final = exp(10 -8 ) [20[. Of course, at each step the 
preservation of the knot topology is checked with the 
method explained before. We show as an example the 
obtained densities of states for the unknot, the trefoil 
and the figure-eight of length L = 100 in fig. HJ 

The obtained distribution of density of energy states 
can now be used to calculate the canonical averages of 
the observables according to standard relations. As al- 
ready mentioned, the observables studied in this paper 
are the internal energy of the canonical system and the 
heat capacity as a function of the temperature. The ex- 
pression of the internal energy is given by: 



W)> = 



(4) 



where m is the number of the closest contacts of polymers 
and Q(E) is the density of energy states computed by 
means of the WL algorithm. The formula of the heat 
capacity has been already provided in Eq. (pQ). 

The computational method has been checked using the 
trivial knot U as a test. We could repeat the previous 
results of [24| in the case of chains with lengths equal to 
12, 30 and 50. 
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FIG. 5: Specific energy and heat capacities for the unknot polymer as functions of the normalized temperature T = The poly- 
mer length can take the values L = 100 (upper triangles), 150 (circles) , 200 (rectangles) , 250 (diamonds), 300 (down-triangles), 
where L denotes the number of steps, (a) Plot of the specific energy (in e- units) of the unknot; (b) Plot of the specific heat 
capacity (in e- units) of the unknot. 



III. RESULTS AND DISCUSSION 

In this section we study the thermal properties of 
the trefoil and figure-eight with lengths in the range of 
100 — 300 segments for the case e > 0. The used algo- 
rithm is quite fast. Depending on the length, the typical 
simulation times vary between 15 minutes to 1 hour on 
a dedicated pc. It can be applied to much more compli- 
cated knots. As an example, we have analysed here also 
the case of the 8i knot. 

The results concerning the trefoil and figure-eight are 
presented in fig. and fig. [3 respectively. For comparison, 
fig. [5] represents the results for U with lengths ranging 
from 100 to 300. It is found that the growth of the spe- 
cific energies (E(/3))/L in fig. EJa) and fig. Ufa) is charac- 
terised by three regions both for the trefoil and the figure- 
eight. At very low temperatures the energy increase is 
practically zero. The reason is that the temperatures are 
too low to allow contacts between the monomers. When 
the energy of the system grows - let's remember that we 
are working in energy units, so that we can talk about 
temperatures and energies interchangeably - the "first 
energy state" can be reached. This consists in the case in 
which at least two monomers are in contact. For instance, 
we expect that when the chain has length L = 100, the 
first energy state occurs in figs. [6] and [71 at some "critical" 
temperature T\ ~ 1, while when L = 200 this temper- 
ature should be given by T\ ~ 2. These values of T\ 
slightly increase with increasing topological complexity. 
In general, this proportionality temperature/length is re- 
spected in figs. [6] and After the first energy state is 
reached, the specific energy grows rapidly as a function of 
the temperature until saturation is reached after a given 



temperature T2 and the energy increase becomes moder- 
ate. The division into three regions can be expected by 
looking at the Boltzmann factor e~^ m£ appearing in the 
expression of the energy given in Eq. (j4j): 

i) In the range T <T±, the factor e~ l3rn£ goes rapidly to 

0, because low temperatures cause few states to be 
occupied in this situation, so we get (E) — >• 0. 

ii) When T ± < T < T 2 , the factor e~^ m£ chan ges sensi- 

tively depending on /? = T _1 . In this regime more 
and more states of the given system are getting 
excited through heat absorption and (E) increases 
with the increase of the temperature T as expected. 

hi) Once T > T2, (E) increases slower due to the fact 
that most states have been already occupied in the 
range T\ < T < T2. Just few states are left that can 
be excited in this situation. One could notice that 
the factor e~^ m£ in Eq. (j4j) increases slowly at high 
temperatures. When T — >> 00, we have e~^ m£ — )► 1. 

It can also be found that in fig. [6fb) and fig. [7(b), the 
heat capacity has a peak in the range of temperature 
T\ < T < T2 . The appearance of this peak together with 
the sigmoidal behavior of the specific energy points out at 
an ongoing two-state phase transition according to [25|. 
The presence of these pseudo-transitions related to the 
lattice geometry have been clarified in ref. [26|, [27| . We 
would also like to show that the heat capacity and en- 
ergy are affected not only by the temperature, but also by 
the polymer length L and its topological conformation. 
First of all, the specific energy of the chain decreases 
with increasing values of L. For instance, the quantity 
(E(f3))/L for L = 100 is larger than that of longer chains 
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FIG. 6: Specific energy and heat capacities for the trefoil as functions of the normalized temperature T = ~. Separated plots 
are provided for trefoils with L = 100 (upper triangles), 150 (circles), 200 (rectangles) , 250 (diamonds), 300 (down-triangles), 
(a) Plot of the specific energy (in e-units) of the trefoil; (b) Plot of the specific heat capacity (in e-units) of the trefoil. 



as it is illustrated in fig. [6] and fig. H This is because 
longer polymers in the presence of short-range repulsive 
interactions prefer states which are diluted comparing 
with the short ones. Probably for this reason, in proteins 
knots are localized in regions with a very small number of 
residues Q. The total internal energy, instead, decreases 
with decreasing polymer lengths. This is due to the fact 
that longer polymers can accommodate a larger number 
of contacts. Let us note that the peak temperatures of 
the heat capacity are shifting with the length L. This is 
expected since we are computing the specific heat capac- 
ity. If we computed the total heat capacity of the system, 
these peaks would occur at the same temperature. 

Concerning the topological effects, we observe that 
both the energy and the heat capacity increase with 
increasing knot complexity. This can be seen by com- 
paring the plots for various knots of the same length in 
figs. [5j [6] and [3 On the other side, the topological effects 
are stronger when polymers are shorter. This is shown in 
Table U in which the energy differences for polymers with 
the same length but different topology are listed. In the 
table dur is the difference between the unknot and the 
trefoil defined by dur = , while dq-j denotes 

the difference between the trefoil and figure-eight. When 
the length ranges from 100 to 300, dur varies from 
0.073 to 0.015 and d TT from 0.032 to 0.004 at T = 10e. 
The heat capacity exhibits a similar behavior. It is 
thus pretty clear that the topological effects become 
less important when the length is increasing. This is 
similar to the so called 'size effect' in nanomaterials (28^ . 
Topological properties will have a lesser influence in the 
case of diluted monomers. We expect that topological 
effects will disappear if the polymer is long enough. 

To get a feeling of the size of the involved quantities, 



TABLE I: Energy differences between the unknot and the 
trefoil as well as those of the trefoil and figure-eight. All knots 
have the same length. The values of the energy difference are 
reported in the case of various lengths L = 100, 200, 300. Here 
the temperature is T = 10s. 
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unknot trefoil fi£ 


^ure- eight 


dur 


d rj r 


100 


0.174 


0.247 


0.279 


0.073 


0.032 


200 


0.170 


0.197 


0.210 


0.027 


0.013 


300 


0.167 


0.182 


0.186 


0.015 


0.004 



we can make the concrete example of the proteins which, 
as mentioned before, can be found in different topological 
conformations like the trefoil and the figure-eight. These 
proteins have a length of about 250 residues [2] and their 
temperature is around 310 K (since this is approximately 
the temperature of the human body). The parameter e 
for a protein is about lOK/mole [29]. Based on those 
data, for one mole, we can roughly estimate the range of 
corresponding temperatures in our units, which is T/e = 
ZIOK/e = 31. This temperature is located in the region 
T > T2. Before concluding this section, we would like 
to stress that the used method for preserving the knot 
topology after a pivot move is very time efficient and 
can be applied to any knot configuration independently 
of its complexity. For example, in fig. [8] we provide the 
data for the specific energy and heat capacity in the case 
of the knot 81. As it is possible to check, also for this 
knot the same conclusions drawn before for /C = T, T 
are valid. In particular, the topological effects become 
less important with increasing polymer length and it is 
confirmed that both the energy and the heat capacity 
increase with increasing knot complexity. 



6 





0.3r 




0.25 - 




0.2- 




0.15- 


LLJ 




0.1 - 




0.05- 








°0 



(a) 



(b) 



150 
200 

250 / 
300 



" 

/ . 

* ■ ♦ T T 



_i i i I i i i_ 



4 6 8 10 
T 





0.16: 




0.14^ 




0.12^ 




0.1 r 


C/L 


08 '- 




0.06: 




0.04^ 




0.02: 




0^ 



— I 1 — I — I — I — I — I — I — I — I — I 1 — I — I — I — I — I — I — I— 








8 10 



FIG. 7: Specific energy and heat capacities for the figure-eight as functions of T= ^. The lengths are L = 
100 (upper triangles), 150 (circles) , 200 (rectangles) , 250 (diamonds), 300 (down-triangles) steps, (a) Plot of the specific energy 
(in e- units) and (b) plot of the specific heat capacity (in e- units) of the figure-eight knot. 
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FIG. 8: Specific energy and heat capacities for AC = 81 as functions of T = ^. The studied lengths are L — 
100 (upper triangles), 150 (circles) , 200 (rectangles) , 250 (diamonds) and 300 (down-triangles) steps, (a) Plot of the specific 
energy (in e- units) and (b) specific heat capacity (in s- units) of the 81 knot. 



IV. CONCLUSIONS 

We have studied the topological effects on the thermal 
properties of several knot configurations with the help 
of Monte Carlo simulations based on the Wang-Landau 
algorithm and the pivot method. The topology of the 
knots, that can change after a pivot move, has been pre- 
served by explicitly checking if during the move some 
parts of the trajectory of the knot have crossed them- 
selves at some points. 

The specific energy and the heat capacity of the trefoil, 
the figure-eight and the 81 knot have been calculated for 



various lengths ranging from 100 to 300. As a test of 
our methodology, these quantities have been derived also 
in the known case of the unknot. In particular, we have 
been able to reproduce the results of [24] . Our investiga- 
tions show that the specific energy and heat capacity for 
a polymer knot whose monomers are subjected to short- 
range repulsive interactions are depending on the length 
and the topological conformation of the polymer knot. 
This could be expected because in the presence of repul- 
sive interactions the polymer will prefer "diluted" con- 
formations in which the monomer density is lower. We 
may thus suggest that the dilution of polymers in a good 
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solvent should be effective in reducing the topological ef- 
fects. Moreover, for the growth of the specific energy 
with respect to the temperature it is possible to distin- 
guish three different regimes that have been explained 
in details in Section II. The zero energy regime appears 
at extremely low temperatures and it is characterised by 
the fact that the energy of the thermal fluctuations is be- 
low the threshold of the repulsive energy barrier, so that 
the monomers have not enough energy to get near. After 
this threshold has been passed, the number of contacts 
between the monomers starts to increase rapidly until at 
the end a saturation phenomenon occurs and the energy 
of the system increases only slowly with increasing tem- 
peratures. The plot of the specific heat capacity shows 
a pseudo-phase transition which is related to the lattice 
geometry, as explained in refs. [26] and [27j. 

Our simulations can be generalised under several 
aspects. First of all, we limited ourselves to simple 
short range interactions, but there are no obstacles 
to extend our procedure to describe more realistic 
polymer systems. Moreover, we have limited ourselves 
to pivot moves involving a short number of segments. 
No difference between the case N = 4 and N = 5 
have been detected, so that all displayed diagrams are 



related only to the case N = 5. The independence of 
N up to N = 10 has been confirmed by preliminary 
calculations performed using another method to preserve 
the topology based on knot invariants. We hope to be 
able to report on these new developments very soon (30| • 
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